Bernoulli potential in type-I and weak type-II superconductors: Surface charge 
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The electrostatic potential close to the surface of superconductors in the Meissner state is discussed. 
We show that beside the Bernoulli potential, the quasiparticle screening, and the thermodynamic 
contribution due to Rickayzen, there is a non-local contribution which is large for both type-I and 
weak type-II superconductors. 



I. HISTORICAL REVIEW 

The electrostatic potential in superconductors is 
known already for seven decades. Within the pre-London 
approaches based on the ideal charged liquid, the elec- 
trostatic potential was assumed to balance the Lorentz 
and the inertial forces acting on diamagnetic currents. 1 
The current j = ertv flows along the surface and its am- 
plitude falls off exponentially from the surface into the 
bulk. Similarly the screened magnetic field decays. The 
balancing electrostatic potential was found to be of the 
Bernoulli type, e<p = — |rm; 2 , therefore it falls off on the 
scale of half the London penetration depth. From the 
Poisson equation, 



-e V p = p, 



(1) 



one finds that this corresponds to a charge accumulated 
in the layer penetrated by the magnetic field. To main- 
tain the charge neutrality, the accumulated bulk charge 
has to be accompanied by the opposite surface charge. 
This picture has been confirmed by the London theory. 2 



A. Advanced London-type approaches 

The electrostatic potential equals the Bernoulli poten- 
tial only at zero temperature when all electrons are in 
the superconducting condensate. At finite temperatures, 
a part n n of electrons remains in the normal state, while 
the rest n s = n — n n contributes to the supercurrent, 
j = en s v. The electric force acting on normal electrons 
is transferred to the condensate via a mechanism which 
reminds one of the fountain effect in superfluid helium. 
As a consequence, the electrostatic potential reduces to 
e<p = -^|mu 2 . This reduction of the potential derived 
by van Vijfeijken and Staas 3 has a somehow confusing 
name: the quasiparticle screening. It should be noted 
that the electrostatic potential caused by currents in su- 
perconductors is traditionally called the Bernoulli poten- 
tial in spite of the quasiparticle screening and other non- 
Bernoulli contributions found later. 

Jakeman and Pike recovered the result of van Vijfeijken 
and Staas from the static and the classical limit of the 



time-dependent Ginzburg-Landau (GL) theory. 4 Unlike 
previous studies, their approach includes the real screen- 
ing due to spatial distribution of the charge, therefore 
it describes the surface charge as a space charge local- 
ized on the scale of the Thomas-Fermi screening length. 
This picture of the surface charge was accepted for a long 
time as it looks natural and seems to be justified from 
the semi-microscopic theory. 

Later studies of the electrostatic potential did not re- 
address the question of the surface charge but attempted 
to overcome the hydrodynamic picture common to all 
the above mentioned studies. Naturally, they focused on 
a detailed description of the Bernoulli potential on the 
scale of the London penetration depth. 

Adkins and Waldram recovered the Bernoulli potential 
from the BCS theory for the system at zero temperature. 5 
They also indicated that at finite temperatures one 
should expect a contribution which depends on the band 
structure. 

The predicted contribution has been evaluated by 
Rickayzen. 6 He used the thermodynamic approach which 
yields 
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The density dependence of the condensate density can 
be evaluated either from the BCS theory leading to the 
result of Adkins and Waldram or from the two-fluid rela- 
tion n s w n (l - 2j ) . The latter gives a simple formula 
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The first term balances the Lorentz and inertial forces, 
while the second reflects the pairing. We will call terms 
proportional to the thermodynamic correction. 

In the next decade, the interest in the electro- 
static potential in superconductors shifted towards non- 
equilibrium situations, see e.g. Ref. 7. It should be noted 
that the Bernoulli potential introduced within the theory 
of non-equilibrium systems is not the same as the one 
we discuss. Within the non-equilibrium theory the term 
Bernoulli potential is used for the velocity-dependent 
part of the energy of Cooper pairs. Only in equilibrium, 
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the velocity-dependent energy agrees with the electro- 
static potential we discuss here. 

B. GL approaches — local 

New studies of the Bernoulli potential in equilib- 
rium have been stimulated by the question of the vor- 
tex charge. Khomskii et al. 8,9 have used the idea of 
van der Marel 10 to derive a simple estimate of the electro- 
static potential in terms of the superconducting gap. The 
obtained formula corresponds to the thermodynamic cor- 
rection of Rickayzen with the original balance term being 
omitted. While the validity of the approach is restricted 
to temperatures close to T c , it has the advantage, how- 
ever, that it can be applied within the GL theory. Blatter 
et al. 11 have employed this estimate to predict a possible 
electrostatic field above the surface of the superconductor 
with the vortex lattice perpendicular to it. 

An opposite limit represents the treatment of the vor- 
tex charge by LeBlanc. 12 He considers only the balance 
term omitting the thermodynamic correction completely 
As one can see from Rickayzen's formula, this approxi- 
mation is limited to the region of low temperatures. 

C. GL approaches — non-local 

All the above mentioned theories have in common that 
the electrostatic potential is a local function of the clas- 
sical kinetic energy or of the gap. In other words, they 
neglect the gradient of the condensate. The non-local 
corrections within the GL theory have been proposed in 
Rcf. 13. It turned out that the non-local corrections are 
taken into account if the classical kinetic energy in the 
balance equation is replaced by its quantum-mechanical 
counterpart. 

For our next discussion it is important that the non- 
local contributions can be rearranged into a local form. 13 
Indeed, within the GL theory the sum of the kinetic en- 
ergy and the GL potential is zero which allows us to ex- 
press the balance term via the local GL potential. Note 
that this rearrangement is possible only within the non- 
local theory, because gradients of the condensate density 
n s = 2\tp\ 2 , where ip is the GL wave function, provide 
important contributions to the quantum kinetic energy, 
in particular close to the surface. 

In Ref. 13, the quasiparticle screening and the ther- 
modynamic correction have been neglected, the com- 
plete theory has been presented in Rcf. 14. The detailed 
derivation together with estimates of material parame- 
ters and the numerical treatment of the Abrikosov vortex 
lattice in niobium at various temperatures can be found 
in Ref. 15. Again, all contributions including the kinetic 
energy can be expressed in a local relation c/?[|-0| 2 ]. 

The non-local form of the Bernoulli potential requires 
to reconsider the picture of surface charges. It is para- 



doxical that the possibility to rearrange the non-local 
corrections into a local function plays the crucial role. 
As noticed already by Yampolskii et al., 16 if the electro- 
static potential is of the form (^[IV'I 2 ], the GL boundary 
condition, Vt/> = in the direction normal to the surface, 
implies that no surface charge is needed. Indeed, at the 
surface V</? = f^V^ = 0. The electric field thus van- 
ishes at the surface which is a sufficient condition for the 
charge neutrality. 

D. Plan of the paper 

Apparently, there is no surface charge on the scale of 
the Thomas-Fermi screening length. The electric field 
due to the dipole between the surface charge and the 
bulk charge was supposed to balance the Lorentz force, 
therefore one has to ask the question how the balance of 
forces looks if the surface charge is absent. 

The charge profile in the vicinity of the surface has 
been numerically studied in Ref. 17. This study was lim- 
ited to thin slabs and type-II superconductors, its results 
indicate, however, that the charge profile at the surface 
tends to organize into two layers - a narrow surface charge 
and a thicker bulk charge. 

In this paper we discuss the charge and potential pro- 
file at the surface of a semi-infinite superconductor in the 
Meissner state. In the next section we shortly describe 
the surface charge within the GL theory. In section III 
we discuss the limit of a weak magnetic field, where all 
results can be derived analytically. We will show that the 
Bernoulli potential at the surface differs from the value 
predicted by Rickayzen. The correction has the form 
of a multiplicative factor that approaches unity for the 
extreme type-II superconductor. Numerical results on 
strong magnetic fields are presented in section IV reveal- 
ing a nonlinear mechanism which, however, is found to 
be very small. The conclusions and outlook are presented 
in section V. 



II. SURFACE CHARGE WITHIN THE 
GINZBURG-LANDAU TYPE THEORY 

A typical profile with the charge distribution close to 
the planar surface and homogeneous external magnetic 
fields is shown in Fig. 1. The profile has been obtained 
by a numerical solution of the GL equation 

ft 2 V 2 e* 2 A 2 jTl ( e V 

2m* ^ + 2m* ^ 2n ^ ^ _ 2^2 J 

(4) 

where n is the density of pairable electrons, 7 is the linear 
coefficient of the specific heat, m* = 2m and e* = 2e are 
the mass and charge of the Cooper pair, and t — T/T c is 
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the temperature on the scale of the critical temperature. 
In the assumed geometry one can use the London gauge 
in which the GL wave function is real. 
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FIG. 1. The charge profile close to the surface for 4 dif- 
ferent magnetic fields. The magnetic field is given in di- 
mensionless units proportional to the magnetic flux quanta, 
Bo = <3?o/(27rAo) = efj,onh/2m. The coordinate x is measured 
in units of the London penetration depth Ao at zero tem- 
perature. The assumed temperature t = T/T c = 0.9 yields 
the penetration depth A = 1.7Ao- The charge is scaled with 
B 2 so that the two lowest magnetic fields B = O.OlBo and 
B = O.lBo give nearly identical curves. The magnetic field 
B = 0.34Bo is close to the critical value. The GL parameter 
kql = 1 corresponds to weak type-II superconductors. 



The density dependence of 7 and T c has been estimated 
from the experimental specific heat, theoretical free- 
electron density of states, and the McMillan formula 19 , 
see Appendix in Ref. 15. The Thomas-Fermi screening 

length, A^ F = - ^ 2 e ° , is very small and its contribu- 
tion can be neglected, because the general integral which 
is proportional to exp(— x/Atf) is zero since no surface 
charge forms on this scale. 

In the final step the charge profile is evaluated from 
the Poisson equation (1). 

III. LOW MAGNETIC FIELDS 

In weak magnetic fields the GL wave function is only 
little perturbed from its bulk value ip = tp^ + Sip with 
Sip <C ip 00 ■ I n this limit the charge profile can be solved 
analytically. 

In the Maxwell equation (5) the vector potential is a 
small quantity so that the perturbation of the wave does 
not contribute in lowest order. The vector potential can 
thus be solved taking ip w ip^ = \/ n (l — * 4 )/2- For 
both types of superconductors, this results in a simple 
exponential decay of the vector potential, 

A x w B\e~ x l x . (7) 

on t he sca le of the London penetration depth A = 
Ao/Vl - £ 4 , with Ag = 2m*/(n a e* 2 n). Here B is the 
value of the magnetic field at the surface at x = 0. 



Equation (4) has been proposed by Bardecn as an ex- 
tension to low temperatures. 18 In fact, Bardeen arrived to 
this equation independently adding the quantum kinetic 
energy to the Gorter-Casimir two-fluid model. Close to 
the critical temperature the GL function has a small mag- 
nitude 2ip 2 <C n which allows one to expand the square- 
root l/^/l — ^ip 2 w 1 + -^ip 2 , an d one finds the custom- 
ary GL equation with parameters a — ^T 2 (t 2 — l)/2n 
and 13 = -fT 2 t 2 /n 2 . 

The vector potential A is given by the Maxwell equa- 
tion 

*2 

V x V x A = u —ip 2 A. (5) 

to* 

Equations (4) and (5) are solved together with 7 and T c 
kept constant. 

From the GL wave function ip we evaluate the Bernoulli 
potential 15 

e V A| F V 2 ^ = -fl ( 1 - * I i\? 



A. GL wave function 

The linear approximation of the GL equation (4) reads 
?i 2 V 2 7 T 2 t 2 / 2 (2 A" 1 ,2 Sl e* 2 A 2 , 

(8) 

which can be expressed in terms of the GL coherence 
length <f = 2fi 2 nt A /{m* 1 T 2 {\ - i 4 )) as 

V 2 # - Aty = ^V'ooA 2 J B 2 c- 2 -/\ (9) 

The solution of (9) is composed from the particular 
integral with the decay on the scale of the London pen- 
etration depth and the general integral decaying on the 
scale of the GL coherence length, 

Sip = tJjxc- 2x/x + ^e _V5x/€ . (10) 

The general integral of (9) can also include a growing 
term cx exp(v / 2x/£). This term is excluded, because the 
perturbation asymptotically vanishes in the bulk, <5-0 — > 
for x — > 00. 



3 



From the GL boundary condition, V6ip — at x = 0, 
we find that the two amplitudes are linked by 



iPz = i>\, 

K 



(11) 



where k — A/£. We note that within Bardeen's exten- 
sion of the GL theory, the ratio A/£ depends on the tem- 
perature, k = K GL /t 2 , where k gl = ^^y^- is the 

ordinary GL parameter defined at T c . The diverging ef- 
fective GL parameter at low temperatures shows that 
within Bardeen's approximation any system behaves as a 
type-II superconductor at sufficiently low temperatures. 
This feature contrasts with the linearly decreasing effec- 
tive GL parameter one finds close to the upper critical 
magnetic field. 15 . Of course, the strong temperature de- 
pendence of k is rather an artifact of the model. Since 
k is the most important quantity of the GL theory, one 
can simply take A and k as input parameters. 

The amplitude of the particular integral is obtained 
from the GL equation (9) as 
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The GL wave function from (10-12) results in 



8ij) = 



e* 2 X 4 B 2 



-2x/X _ 



V2 



2h 2 



(12) 



(13) 



B. Bernoulli potential 

Finally, we take the linear approximation of the 
Bernoulli potential (6) in Sip. Since the asymptotic value 
of the potential (p^ = (^[l^ool 2 ] is not essential, we focus 
on the perturbation caused by the magnetic field 

eSip - A TF V 2 e<V = CSip . (14) 

In terms of A and k the coefficient C reads 
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A general integral which decays into the bulk is of the 
form 

(16) 

From (14) one finds the amplitudes 



The amplitude of the general integral is given by the 
charge neutrality, which requires V(5<p = 0at.T = 0. 
From (16) follows 

eSipTF = -2^ZCa^a - \/2^ZQVe- (19) 

The relative amplitudes proportional to Atf / A and 
Atf/C arc not sufficiently small to be neglected. The 
charge density is given by the second derivatives so that 
one obtains relative amplitudes of the charge density at 
the surface proportional to A/Atf or to £/Atf- The sur- 
face charge is an integral, i.e. the contribution in question 
is proportional to the amplitude at the surface multiplied 
with the Thomas-Fermi length. According to condition 
(19) all charges can have comparable values. 

To prove that the screening on the scale of the Thomas- 
Fermi length Atf and the related surface charge can be 
neglected, one has to employ the GL boundary condition 
\7Stp — at x — 0. Using (11) one obtains 
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so that the general integral is proportional to the cube of 
the small ratio Atf / A. Accordingly, the charge due to the 
general integral is smaller by a factor A TF /A 2 than the 
charge created by the particular integral. In the following 
we neglect this contribution for simplicity of notation. 
Within the same level of accuracy we take C\ = C and 
Q = C. 

Now we are ready to evaluate the electrostatic poten- 
tial. Using (11) and (12) we obtain from (16) 
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This linearized relation extends Rickayzen's formula to 
type-I and weak type-II superconductors. 

Let us first take a look at the extreme type-II su- 
perconductor, k 3> 1, for which Rickayzen's formula 
is recovered. For the assumed real GL wave function, 
the current is proportional to the vector potential j = 
^-ibz^A = — n s A. With the velocity defined via the 
current, j = e*ip 2 <1 v = en s v, one finds from (7) the rela- 
tion between the magnetic pressure and the kinetic en- 
ergy Q~ 2x/X ^^i = ^\ mv ' 1 - Equation (21) yields then 
the Bernoulli potential as 



We have used n s = n(l — t A ) and n n — nt 4 to express 
the temperature dependence in terms of the condensate 
fraction. 

One can see that Rickayzen's formula (3) holds for ex- 
treme type-II superconductors, k — > oo. In this limit the 
general integral has a vanishingly small amplitude and 
the factor 1/(1 + — ) goes to unity. For finite k the gen- 
eral integral deforms the profile of the potential and the 
particular integral has a reduced amplitude. 

Within Bardeen's extension of the GL theory, Rick- 
ayzen's formula is always recovered at low temperatures 
as k diverges with 1/t 2 . Since Rickayzen's formula ap- 
proaches the plain Bernoulli potential for t — > 0, the pre- 
sented derivation also reproduces London's result. 



C. Behavior for finite k 

As already shown, for the extreme type-II supercon- 
ductor, k — ► oo , the general integral vanishes as l//t and 
far from the surface it decays faster than the particu- 
lar integral. The Bernoulli potential thus extends on the 
scale A/2. 

For the finite GL parameter in the region n > \[2 the 
overall picture is rather similar. The general integral is 
appreciable only close to the surface so that the scale 
on which the Bernoulli potential extends is still A/2. Be- 
sides, the amplitude of the particular integral is enhanced 
by the factor 1/(1 - 2/k 2 ). 

For the region k < >/2, the role of the general and par- 
ticular integrals are reversed. First of all, the amplitude 
of the general integral is larger than the amplitude of the 
particular integral. Since the general integral also decays 
more slowly in this case, the Bernoulli potential is dom- 
inated by the general integral and extends on the scale 
£/\/2- The reversed role is also supported by the fact 
that the signs of both parts are opposite as compared to 
the case k > \[2. This sign reversal appears due to the 
enhancement factor 1/(1 — 2/k 2 ). 

Apparently, both components of S(p diverge as k ap- 
proaches \pl. Their sum remains regular, however. For 
this particular case the scale of both contributions is iden- 
tical, namely A/2. It is advantageous to use expression 
(21) to obtain the asymptotic form of the Bernoulli po- 
tential 
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The Bernoulli potentials for different values of the mag- 
netic field are shown in Fig. 2. The limit of low magnetic 
fields agrees with the analytical formula (22), of course. 
One can see that for the two lowest fields B = 0.01 and 
B = 0.1, the potential keeps its profile while its mag- 
nitude scales with B 2 . For magnetic fields close to the 



critical value B c2 w 0.34, the linear theory does not ap- 
ply- 
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FIG. 2. Bernoulli potential as a function of the distance 
from the surface for 4 different magnetic fields. The parame- 
ters are the same as in Fig. 1. 



D. Charge profile 

The charge density results from the Bernoulli potential 
(21) and the Poisson equation (1) as 
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Note that the relative amplitude of the two contributions 
is reciprocal to the relative amplitude of the potentials. 
Integrating from zero to infinity one can check that the 
total charge is zero, as it is required by the charge neu- 
trality. 

The resulting charge density is small compared to the 
total density of electrons en. The ratio of the amplitude 
of (24) to the total electron density can be written as 
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Even for upper estimates of the magnetic field by the 
critical value, B w B c , the factor (25) is extremely small. 
For niobium B 2 / fi Q = ± 7 T C 2 = 3.2 x 10 4 J/m 3 , while 
mc 2 rts 10~ 13 J. This suggests that the maximum charge 
density will be of the order of 10 18 electrons per cubic 
meter. This is by ten orders of magnitude smaller than 
the electron density of niobium n = 2.2 x 10 28 m~ 3 . This 
small value of the accumulated or depleted charge justi- 
fies to ignore the charge profile in the Maxwell equation 
(5) and the change of material parameters like 7, n and 
T c in the GL equation (4). 
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Similarly to the potential, the charge profile reflects 
two scales, A/2 and £/\/2. Before we analyze its proper- 
ties for various values of the GL parameter k, we want 
to discuss a few general features. First, the surface value 
of the charge, obtained for x = in (24), is always neg- 
ative, i.e. the charge carriers are depleted. We call the 
layer of the negative charge the surface charge. Second, 
in the bulk sufficiently far from the surface the charge 
is positive, i.e., the charge carriers are accumulated. We 
call the region of positive charge the bulk charge. Third, 
the width w of the surface is given by p(w) = 0. From 
(24) follows 
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In the extreme type-II superconductor, k — > oo, the 
surface charge is formed by the contribution on the scale 
of the GL coherence length £. The width of the surface 
goes to a value w — > In . The bulk charge extends 
on the scale of the London penetration depth as it is 
known from the classical picture. 

For the limiting case k — > y/2 the width of the surface 
charge is A/ 2. The charge density has the profile 
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Note the different sign inside the last parentheses com- 
pared to the Bernoulli potential (23). Since A/2 = £/\/2 
one cannot associate the surface and the bulk charge to 
the GL coherence length £ or to A. 

For k < \/2, the surface charge is formed by the contri- 
bution on the scale of the London penetration depth. For 
k <C V2, the width of the surface charge is w — ► ^ In . 
The bulk charge extends on the scale of the GL coherence 
length £ as one can see from (24). 



IV. STRONG MAGNETIC FIELD 

In strong magnetic fields the GL wave function is sup- 
pressed in the vicinity of the surface which has to be ac- 
counted for in the Maxwell equation. In this case one has 
to face the fact that the system of equations (4) and (5) 
is nonlinear. We present only a few numerical solutions 
to point out some features typical for strong magnetic 
fields. 

The nonlinear effects in strong magnetic fields follow 
from a deviation of the wave function from its bulk value. 
In the dimensionless representation used in Fig. 3 the 



bulk value is i/'oo = \Jn/2y/\ — t 4 ss 0.6. For the lowest 
field B = 0.01, the deviation is nearly zero on the scale 
of the graph. Also for B — 0.1 the approximation \8tp\ <C 
V'oo is justified. For higher fields the deviation becomes 
appreciable. 
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FIG. 3. Th e deviation o f the wave function from its bulk 
value Voo = ^n(l-t 4 )/2 ^ 0.6 as a function of the distance 
from the surface. The parameters are the same as in Fig. 1. 

The suppression of the wave function at the surface 
shown in Fig. 3 affects the Bernoulli potential by two 
different mechanisms. First, the screening current is re- 
duced at the surface so that the magnetic field penetrates 
deeper into the superconductor. This feature can be seen 
in Fig. 2, however, it is better visible in the profile of the 
charge density shown in Fig. 1, as a widening of the sur- 
face layer. 

Second, the thermodynamic contributions to the 
Bernoulli potential arc nonlinear in the wave function, 
see (6). In contrast, the linearized approximation of the 
Bernoulli potential (21) or (22) does not depend on the 
density derivative of the linear coefficient of the specific 
heat, dj/dn, included in the derivative of the condensa- 
tion energy \^T 2 . Beyond the linearized approximation, 
this derivative contributes as demonstrated in Fig. 4 for 
parameters of niobium. 
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FIG. 4. The charge density which is proportional to the 
density dependence of the linear coefficient of the specific heat 
versus distance from the surface. We use the value for nio- 
bium, d In 7/9 Inn = 0.42, which is close to the parabolic 
band approximation 7 par oc n 1 / 3 , i.e., dln7 par /<91nn = 1/3. 
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Comparing the contribution of dj/dn with the total 
charge density shown in Fig. 1, one can see that this par- 
ticular nonlinear mechanism is rather small giving less 
than 3% of the total charge density. This small value 
justifies the approximations 8 ' 9,11 using only the term pro- 
portional to dT c /dn. 

V. CONCLUSIONS 

In conclusion, we have shown that the Bernoulli 
potential in superconductors can be discussed by the 
charge build up in the region between the surface and 
max[^/y/2, A/2]. The thinner depleted charge region 
carries the surface charge. In contrast to former the- 
ories which assumed that the surface charge is local- 
ized on the scale of Thomas- Fermi screening length Atf , 
we find that the surface charge extends over a range 
L = min[t;/^/2, A/2]. In fact, there is a nonzero con- 
tribution on the scale of Atf, but of a negligibly small 
amplitude proportional to A TF /L 2 . 

For extreme type-II superconductors we have con- 
firmed the picture known already from London with 
the thermodynamic corrections by Rickayzen. The bulk 
charge extends on the scale of the London penetration 
depth A/2. The surface charge is localized on the scale 
of the GL coherence length £/\/2 which is negligible in 
this limit anyway. This is in agreement with the classical 
treatment which neglects all gradient contributions. 

In type-I or weak type-II superconductors, i.e for the 
GL parameter k < \/2, one finds the opposite situation. 
The surface charge is localized on the scale of the London 
penetration depth A while the bulk extends on the scale 
of the GL coherence length £/\/2. In this case the gra- 
dient contributions are dominant and local London-type 
theories naturally fail. 

In strong magnetic fields the Bernoulli potential be- 
comes a non-quadratic function of the magnetic field. 
This feature parallels the nonlinear susceptibility as it is 
dominated by the enhanced penetration of strong mag- 
netic fields into superconductors. There is an additional 
thermodynamic correction absent in low magnetic fields, 
however, for realistic material parameters (like for nio- 
bium) it accounts for less than 3% of the charge density. 
Thus the approximations used in 8,9 ' 11 are justified. 

Finally we remind that the potential <f> discussed here 
does not include the surface dipole which affects the po- 
tential seen outside the sample. We will discuss this con- 
tribution in a forthcoming paper. 
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